#include <iostream>
#include <fstream>
#include "FDmethod.h"
#include <cmath>
#include <algorithm>
#include <vector>
#include <jsoncpp/json/json.h>


class Fun1 : public Function //给定函数
{
public:
    double operator()(double _x, double _y)
    {
        return std::exp(_y+std::sin(_x));
    }
    double diff_x(double _x, double _y) 
    {
        return std::cos(_x)*std::exp(_y+std::sin(_x));  
    }
    double diff_y(double _x, double _y) 
    {
        return std::exp(_y+std::sin(_x));  
    }
    double Laplace(double _x, double _y)
    {
        return std::sin(_x)*std::exp(_y+std::sin(_x)) - std::cos(_x)*std::cos(_x)*std::exp(_y+std::sin(_x)) - std::exp(_y+std::sin(_x));
    }
};

class Fun2 : public Function //自定义函数
{
public:
    double operator()(double _x, double _y)
    {
        return std::sin(_x*_y);
    }
    double diff_x(double _x, double _y) 
    {
        return std::cos(_x*_y)*_y;  
    }
    double diff_y(double _x, double _y) 
    {
        return std::cos(_x*_y)*_x;  
    }
    double Laplace(double _x, double _y)
    {
        return std::sin(_x*_y)*_y*_y + std::sin(_x*_y)*_x*_x;
    }
};

class Fun3 : public Function //自定义函数
{
public:
    double operator()(double _x, double _y)
    {
        return std::exp(-_x*_y);
    }
    double diff_x(double _x, double _y) 
    {
        return std::exp(-_x*_y)*(-_y);  
    }
    double diff_y(double _x, double _y) 
    {
        return std::exp(-_x*_y)*(-_x);  
    }
    double Laplace(double _x, double _y)
    {
        return -std::exp(-_x*_y)*_x*_x-std::exp(-_x*_y)*_y*_y;
    }
};

int main()
{
    Json::Reader reader;
    Json::Value root;
    std::ofstream fout;
 
	std::ifstream in("main.json", std::ios::binary);
	if (!in.is_open())
	{
	  std::cout << "Error: cannot open file." << std::endl; //报错：无法打开文件
	  return -1;
	}

    std::vector<double> h;
    int function_number, boundary_number;

    if (reader.parse(in, root))
	{
        function_number = root["function_number"].asInt();
        boundary_number = root["boundary_number"].asInt();
        for (int i = 0; i < root["n"].size(); i++)
	    {
	        int temp_n = root["n"][i].asInt();
	        h.push_back(1.0/temp_n);
	    }
	}

    fout.open("results.txt");
    fout << "Results of Function" << function_number << ", condition " << boundary_number << ": " << std::endl;

    for (int i = 0; i < h.size(); i++)
	{
            if (function_number == 1)
            {
                Fun1 F;
                FD_domain1 FD(F, h[i], boundary_number);
                FD.solve();
                std::vector<double> error = FD.norm_error();
                fout << "n = " << static_cast<int>(1.0 / h[i]) << std::endl;
                fout << "error norms(1, 2 and infinity) are " << error[0] << ", " << error[1] << ", " << error[2] << std::endl;
                fout << "results(x-positive direction, y-positive direction): " << std::endl;
                for (int j = 0; j <= static_cast<int>(1.0 / h[i]); j++)
	            {
                    for (int k = 0; k <= static_cast<int>(1.0 / h[i]); k++)
	                {
			  if (k == 0)
			    {
			      fout << "[" ;
			    }
			  fout << FD.result(k*h[i], j*h[i]);
			  if (k == static_cast<int>(1.0 / h[i]))
			    {
			      fout << "], " << std::endl;
			    }
			  else
			    {
			      fout << ", ";
			    }
			}
		    }
                fout << "error at (0.125,0.125),(0.875,0.875),(0.125,0.875),(0.875,0.125): " << std::endl;
                fout << std::abs(FD.result(0.125,0.125) - F(0.125,0.125)) << ", ";
                fout << std::abs(FD.result(0.875,0.875) - F(0.875,0.875)) << ", ";
                fout << std::abs(FD.result(0.125,0.875) - F(0.125,0.875)) << ", ";
                fout << std::abs(FD.result(0.875,0.125) - F(0.875,0.125)) << std::endl;
	    }
            if (function_number == 2)
            {
                Fun2 F;
                FD_domain1 FD(F, h[i], boundary_number);
                FD.solve();
                std::vector<double> error = FD.norm_error();
                fout << "n = " << static_cast<int>(1.0 / h[i]) << std::endl;
                fout << "error norms(1, 2 and infinity) are " << error[0] << ", " << error[1] << ", " << error[2] << std::endl;
                fout << "results(x-positive direction, y-positive direction): " << std::endl;
                for (int j = 0; j <= static_cast<int>(1.0 / h[i]); j++)
	            {
                    for (int k = 0; k <= static_cast<int>(1.0 / h[i]); k++)
	                {
			  if (k == 0)
			    {
			      fout << "[" ;
			    }
			  fout << FD.result(k*h[i], j*h[i]);
			  if (k == static_cast<int>(1.0 / h[i]))
			    {
			      fout << "], " << std::endl;
			    }
			  else
			    {
			      fout << ", ";
			    }
			}
		    }
                fout << "error at (0.125,0.125),(0.875,0.875),(0.125,0.875),(0.875,0.125): " << std::endl;
                fout << std::abs(FD.result(0.125,0.125) - F(0.125,0.125)) << ", ";
                fout << std::abs(FD.result(0.875,0.875) - F(0.875,0.875)) << ", ";
                fout << std::abs(FD.result(0.125,0.875) - F(0.125,0.875)) << ", ";
                fout << std::abs(FD.result(0.875,0.125) - F(0.875,0.125)) << std::endl;
            }
            if (function_number == 3)
            {
                Fun3 F;
                FD_domain1 FD(F, h[i], boundary_number);
                FD.solve();
                std::vector<double> error = FD.norm_error();
                fout << "n = " << static_cast<int>(1.0 / h[i]) << std::endl;
                fout << "error norms(1, 2 and infinity) are " << error[0] << ", " << error[1] << ", " << error[2] << std::endl;
                fout << "results(x-positive direction, y-positive direction): " << std::endl;
                for (int j = 0; j <= static_cast<int>(1.0 / h[i]); j++)
	            {
                    for (int k = 0; k <= static_cast<int>(1.0 / h[i]); k++)
	                {
			  if (k == 0)
			    {
			      fout << "[" ;
			    }
			  fout << FD.result(k*h[i], j*h[i]);
			  if (k == static_cast<int>(1.0 / h[i]))
			    {
			      fout << "], " << std::endl;
			    }
			  else
			    {
			      fout << ", ";
			    }
			}
		    }
                fout << "error at (0.125,0.125),(0.875,0.875),(0.125,0.875),(0.875,0.125): " << std::endl;
                fout << std::abs(FD.result(0.125,0.125) - F(0.125,0.125)) << ", ";
                fout << std::abs(FD.result(0.875,0.875) - F(0.875,0.875)) << ", ";
                fout << std::abs(FD.result(0.125,0.875) - F(0.125,0.875)) << ", ";
                fout << std::abs(FD.result(0.875,0.125) - F(0.875,0.125)) << std::endl;
	    }
	}

    fout.close();

    return 0;
}
